Photon transport simulation using the Monte Carlo method in Rust

Metadata
aliases: []
shorthands: {}
created: 2022-07-02 13:26:50
modified: 2022-07-02 14:30:00

This document contains the description and brief documentation of a photon transport simulation program, which uses the Monte Carlo method to simulate the path of particles and is written in Rust.

The program is open source and the source code is available on GitHub.

WASM build scripts are from the eframe_template repository.

User manual

Simulates the transport of monoenergetic gamma photons emitted from a point-like emitter through a scintillation detector and outputs the spectrum made by the photons. The detector is a cylinder shaped Sodium Iodide (NaI) detector centered in the origin.

The program can be run as a native app or available on the internet as a webapp, at this link.

Screenshot of the web build

After opening, the settings panel is visible on the right side of the display. Its use is self-explanatory. There is also a button here, with which we can start the simulation with the set values, and then we can also stop it with this button.

During/after running the simulation, you can also see extra information and statistics about the current run. After stopping the simulation, the results can be saved to a text file with the Export to JSON button (this function is not available on the web interface due to limitations).

The structure of the program

The main and start functions that start the app are located in the files main.rs and lib.rs. The main function is used by the native version, and the start function is used to start the web version.

In Monte Carlo simulations, the computer spends a significant part of the calculations generating random numbers, which is why it is worth using the most efficient random number generating algorithms. Accordingly, I used the WyRand algorithm, which the program imports from an external source in the rand_gen.rs file. The native and web builds had to use separate packages for compatibility, but both use the same algorithm (see fastrand and nanorand).

The Vector structure implemented in the vec3.rs file enables easy handling of 3D vectors. In addition to the basic vector operations, there is also a function that generates a random vector with isotropic direction distribution (random_isotropic_normed), and also one that rotates the vector by a given angle in a random direction (rotate_random_by_angle).

The photon.rs file implements the structure representing photons. Here you can also find some global variables describing the detector, as well as the function that calculates the effect cross-sections. The photons can be simulated along their entire path with the simulate subfunction. To sample the Compton scattering, we use the Kahn method, which is implemented in the sample_kahn_method function.

The photon_emitter.rs file contains code fragments for launching photons in a directed cone to save resources.

The code that organizes the logic of the program and draws the user interface can be found in the app.rs file. Here, the MyApp structure contains all the data of the program. The associated update function draws the user interface, for which the program uses the egui/eframe framework.

When starting the simulation, the function start_simulation is called, which starts the execution of calculations on several threads. Due to the use of multiple threads, the program stores a lot of data in so-called "atomic" variables, thereby avoiding race conditions. In the case of the web version, it is not possible to use multiple threads, so the simulation runs inside the update function on the main thread (this is why the web version is quite slow).

Compilation

The program can be compiled with the standard command in Rust:

cargo build --release

And the same is true for running it:

cargo run --release

Example spectra

Let's examine two spectra that were outputted by the program. Let's call them spectrum and . They used the following settings:

Parameter spectrum spectrum
Position of photon emitter (3;-3;2) (4;0;0)
Energy of photons [] 661,7 1332,5
Radius of detector [] 2,5 3,0
Height of detector [] 3,0 5,0
Density of detector [] 3,67 3,67
FWHM [] 6,0 8,0

spectrum

The spectrum is shown in the following graph:

A spectrum|500

The spectrum shows several characteristic features, which are characteristic of the spectra recorded by scintillation detectors from gamma photons. One of these is the well-known Compton-edge, which is approx. It is located at the end of the plateau from to . Naturally, Compton scattering is responsible for this. The other well-known shape that can be seen is the total-energy peak, which is located at . This is created after photons are completely absorbed through the photoeffect.

The spectrum was created from the simulation of photons, which caused hits.

spectrum

The spectrum is shown in the following graph:

B spectrum|500

The Compton edge is also visible on this spectrum, which now became due to the higher energy of the photons. The full-energy peak can also be found here, which is also located at the initial energy of the photons. In this case, this meant .

In contrast to spectrum , the signs of pair production appeared on spectrum , in the form of two smaller peaks, which are found at and energies. They can occur because the energy of the initial photon is greater than , which allows pair production. The peak with lower energy belongs to the case when both photons escape after annihilation, and the peak with higher energy to when only one photon escapes.

The spectrum was created from the simulation of photons, which caused hits.

Theoretical investigation of Compton edges

According to theoretical calculations, in the case of the spectrum, the Compton edge should be at , while in the case of it should be at . These values differ minimally (by less than 1%) from those measured on the spectra obtained from the simulation, which is good news for us: it shows that our simulation program works effectively in terms of approximating reality.

Between the Compton edge and the total energy peak

If we look at the energy range between the Compton edge and the full energy peak, here again we can see many bumps in both spectra. These impacts are created by photons Compton-scattered multiple times, which eventually escaped from the detector volume.

Change of efficiencies while moving the emitter

With the settings corresponding to the spectrum, we can see the efficiencies recorded during the uniform movement of the photon source in the following graph:

Efficiencies while moving the photon emitter|500

The source was moved between the points with coordinates (1.0; 3.5; 2.0) and (-4.0; -1.5; 2.0), evenly, by recording spectrum. In each measurement, we simulated at least photons, but usually a multiple of this.

During the movement, the source first approached the detector and then moved away from it in the second half. Accordingly, the total efficiency was higher at the closer points and lower at the more distant ones. This is due to the fact that in the closer states the detector covers a larger part of the source's "angle of view", so more photons collide with it, giving it more energy.

However, the intrinsic efficiency decreased as a result of the approach. This could be due to the fact that during the approach, the photons passing through the detector traveled a shorter and shorter distance on average within the detector, since the photons starting from closer and arriving at the edge of the detector base only pass through the corner of the detector, which means a smaller distance, compared to as if they were going along almost parallel to it, which happens with relatively more photons when the source is further away. So, due to the reduction of the average path length traveled within the detector, the individual photons can emit less energy when the source is closer, thus reducing the intrinsic efficiency. Moreover, with the increase in the number of experiments, a symmetrical figure centered on experiments 4-5 was formed due to the rotational symmetry of the detector.

Change of efficiencies while increasing energy

With parameters corresponding to the spectrum, the efficiencies recorded during the uniform increase of the energy of the photon source can be seen in the graph below:

Efficiencies while increasing the energy of the photons|500

The energy was increased between and , in increments of , by performing 11 measurements. Each spectrum was generated by launching at least photons.

It cannot be clearly seen with the naked eye, but here the intrinsic and total efficiencies change in the same way, their ratio remains essentially constant depending on the energy. This is due to the fact that the location of the source did not change between experiments.

Ratio of intrinsic and total efficiencies|500

It is also worth noting here that efficiencies decrease with increasing energy. This is due to the fact that at higher energies the sum of the cross-sections is smaller, due to which the probability of interaction also decreases. As a result, the photons have less chance to release their energy.

The different cross-sections as a function of photon energy|500

Based on the figure, it can be said that the evolution of the efficiencies according to energy follows quite well the change of the sum of the efficiency cross-sections. The decrease is caused by a strong decrease in the sum of the cross-sections for Compton scattering and the photoeffect, in comparison, pair production only increases above the corresponding threshold (), but its increase is too slow to balance the other two effects.

What is important to note, however, is that the sum of the cross sections is not completely decisive, because the nature of each interaction is also important: for example, the photoeffect leads to complete absorption in all cases, but Compton scattering and pair creation allow further photons to escape. For these reasons, the cross-section of the photoeffect plays a greater role compared to the others.